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The relation between the shape of the force driving a turbulent flow and the upper bound 
on the dimensionless dissipation factor (3 is presented. We are interested in non-trivial 
(more than two wave numbers) forcing functions in a three dimensional domain periodic 
in all directions. A comparative analysis between results given by the optimization prob- 
lem and the results of Direct Numerical Simulations is performed. We report that the 
bound on the dissipation factor in the case of infinite Reynolds numbers have the same 
qualitative behavior as for the dissipation factor at finite Reynolds number. As predicted 
by the analysis, the dissipation factor depends strongly on the force shape. However, the 
optimization problem does not predict accurately the quantitative behavior. We com- 
plete our study by analyzing the mean flow profile in relation to the Stokes fiow profile 
and the optimal multiplier profile shape for different force-shapes. We observe that in our 
3D-periodic domain, the mean velocity profile and the Stokes fiow profile reproduce all 
the characteristic features of the force-shape. The optimal multiplier proves to be linked 
to the intensity of the wave numbers of the forcing function. 



1. Introduction 



Richardson 


(1922 


), 


Taylor 



( 1938 1 and Kolmogorov ( 1941 1 developed a description of turbulence based on the concept 
of an energy cascade. In this defining description of the turbulence phenomenon, kinetic 
energy is transferred at a constant rate from larger unstable eddies to smaller eddies. The 
cascade lasts until the eddy motion becomes stable and viscosity can effectively dissipate 
the kinetic energy. The rate of dissipation of kinetic energy s is defined by: 



£ = 2iy{SijSi. 



(1.1) 



where Sij is the rate of strain tensor and v the kinematic viscosity coefficient. Analyses 
of the Navier-Stokes equations show that bounds on the kinetic energy dissipation rate 
can be derived directly from these governing equations without additional hypothesis or 
approximations for statistically stationary incompressible ffows. During the 1990s, rigor- 
ous bounds for different types of boundary driven flows were derived in the asymptotic 



case of an infinite Reynolds number Re (Doering & Constantin 1992 Marchioro 1994 
VV^[T997| : 

Re^Ul/v, (1.2) 
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where U — \J (u^) is the steady-state root mean square value of the total velocity field and 
I a characteristic length scale of the flow. Variational methods for optimizing the values 



of the bounds were later introduced by Doering & Constantin (1994) and Constantin h 



Doering ( 1995 1 . Since then, these methods have been at the center of the quantitative 



evaluation of the upper bounds on the dissipation rate in boundary-driven turbulent flows 
( Nicodemus et aL]|1998 Kerswell]|1998 ). The first rigorous limits on bulk dissipation for 
body-force-driven turbulence in a fully periodic domain were derived by IChildress t t al.\ 
( [2001| ) and Doering & Foias ( 2002 1 . Their result indicates that for the three-dimensional 
Navier-Stokes equations, the energy dissipation rate e satisfies: 



(1.3) 



where I is the longest characteristic length scale in the body-force function and c\ and C2 
are coefficients that depend only on the functional shape (defined more precisely in the 
next section) of the body-force. Dividing equation 1.3 by jl yields: 



/? < -f C2, 

tie 



(1.4) 



where /3 = el/U^ is the dimensionless dissipation factor. The behavior of the dissipation 
factor implied by this relation is in qualitative agreement with theoretical, computational 
and experimental results for homogeneous isotropic turbulence ( Frisch|1995, : Sreenivasan] 



1984 19981. Doering et at. (2003) derived an explicit upper bound for (5 depending 



only on the shape of the external forcing at high Reynolds number. The perspective of 
determining a value for e, via /3, directly from the external forcing function presents 
invaluable advantage since nearly all current turbulent models rely on the prediction of 
£ from a solution of a transport equation or semi-empirical models. 

In this short paper, we present the results of a systematic study of the qualitative 
features of the bounds on the dissipation factor. Few comparisons have been made be- 
tween the bounds obtained through analytical derivations and the results given by Direct 
Numerical Simulations (DNS). Moreover, only canonical large scale forcing such as Kol- 
mogorov flow (Childress et aZ.|2001l or constant shear (Doering et al. 2003) has been 



tested. Using non-trivial forcing functions, we focus, first, on the extent of the (3 depen- 
dence on the body- force-shape. Second, we analyze the mean profile dependence on the 
force profile. 



2. Bound on the dissipation factor 

2.1. Definitions 

We are considering a body- force driven incompressible flow on a minimal domain ( Jimenez 
fc Moin|[i991) periodic in three directions. The dynamics of the flow is governed by the 



Navier-Stokes equations for an incompressible fluid 

— + (u- V)u + Vp=j.V2u + f (2.1) 

V-u = 0, (2.2) 

where u is the velocity vector, p is the pressure and f the driving force. The derivation of 
the bound on the dissipation factor /3 in our case is identical to the one made by [Doeringl 



et al. (2003) for the case of a body-forced plane shear flow, except for the boundary 
conditions. We therefore limit ourselves to recalling the key points of that derivation and 
refer to their work for the detailed description of the solution of the variational problem. 
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The steady force driving the fluid can be written as: 

f (x) = F0(y/Oe,, (2.3) 

where I is the longest length scale in the forcing function, here the characteristic length 
of the domain, and F the amplitude. The dimensionless square integrable (or smoother) 
shape function (f) E i^[0, 1] satisfies homogeneous Neumann boundary conditions with 
zero mean: 

^'{0) = cb'{l), [ ,^(7;)d77 = 0. (2.4) 



The shape function is normalized by: 

1= / m'dv- (2.5) 



1 

2, 





defining an unique amplitude F for a given f . For practical purposes, the dimensionless 
potential $ S H^[0,l], the space of functions with square-integrable first derivatives, 
is introduced. = —0 and <i> satisfies homogeneous Dirichlet boundary conditions, 
$(0) = $(1). 

Next, we introduce a mean zero multiplier function t/j G _ff^[0, 1] not orthogonal to 
(j), {(jiip) 7^ 0, and satisfying homogeneous Neumann boundary conditions ip'{0) = i/j'i^)- 
We define, ^' e i/^[0,l], the derivative of the multiplier function (^ = ip') satisfying 
homogeneous Dirichlet boundary conditions ^'(0) = ^'(1)- The inner product of (j) and ^jj 
is equal to the inner product of <i> and i.e. {^'i') — {(pijj) ^ 0. 

2.2. Variational problem and solution 

For a steady state flow, the mean energy injected by the forcing F{(j>Ux) should be equal 
to the total dissipation e : 

e = iy{\Vu\^)=F{cl)u,). (2.6) 

Another expression for the forcing amplitude F can be obtained by projecting the 
streamwise component of the Navier-Stokes equations onto the multiplier function tp. 
The inner product of the Navier-Stokes equations with il^iij/Vje^ is integrated by parts 
over the volume. Then the long-time average yields: 



(2.7) 



Then equation |2.6| becomes 



An expression for the dimensionless dissipation factor (3 is obtained by dividing by U'^ /I: 

Changing the velocity variables to normalized velocities ue^ + vBy + we^ ~ U^^{uxex + 
e^), so that {u^ + v''^ + w^) — 1, and using the potential <i> and derivative 
multiplier 'i', equation |2.9| is recasted as: 

P- ^^-^"^ 
The upper bound /3h on the dissipation factor is obtained by maximizing the right-hand 
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side of |2.10| over the normalized velocity field, and then minimizing over all multiplier 
functions 



P ^ min max 



(cf>'ii) (vpuu) 



Re 



_i($'ii)(*'u> 



(2.11) 



for any solution of the Navier-Stokes equations. depends explicitly only on the Reynolds 
number and on the shape (4>' = —</)) of the applied force. 



The evaluation of /3b(i?e) going beyond the scope of this paper, we refer to [Doering 
:oT that matte 

l3b{Re) — min 



et al. ( 2003 1 for that matter and just give their result 

2^1/2 

n — 



1 



27 ye[o,i 



sup I ^(y) I +i?e-i(*'2)i/2 



(2.12) 



Finally, these authors solve exactly the extremization problem for the optimal in the 
asymptotic case Re — > oo: 



lim sup /3&(i?e) ^ Pb{co) = 



(2.13) 



The shape function being normalized (i.e. 
further: 

/36(oo) = 



1), we can simplify this bound a little 



(2.14) 



27(1 $ I) 



3. Force-shape dependence of the dissipation factor 

A non-trivial body-force-shape is tested based on the Kolmogorov forcing. A second 
wave number is added with different amplitudes Ak to the wave number k = 1 used in 
the traditional Kolmogorov forcing: 

/(x) = [sin(27rr;) + sin(2^fc7?)] e, = ^,^(77) e„ G [0, 1] , fc > 1 (3.1) 

where. 



F = \ — ^ and 0(77) = J [sin(27rr/) + Afc sin(27rfc?7)] , (3.2) 



1 + 4 



2 

with ^1 = and g 3? for fc ^ 2. In the following, the term representing the traditional 
Kolmogorov forcing in equation 3.1 will often be referred as "primary term" and the term 
of higher wave number as "secondary term" . 

We can now use |2.14| in order to evaluate the bound on /3 in the asymptotic case 
Re — > 00 for the non-trivial body-force-shape define above: 

/3fc(00) = ^ — ; . (3.3) 

^ /o I 2^\/t5 (cos(27r7y) + ^ cos{2nkv)) \ drj 



The integral on the denominator of the right hand side of equation |3.3| can be rather 
painful to solve by hand. First, we consider the simple cases. 

The case of the classical Kolmogorov forcing, fc = 1 {Ai = 0), is straightforward. 
Equation 3.3 becomes: 

1 1 77^ 

Pb{oo) = = (3.4) 

^ Jo I cos{2n7i) I d,7 V54 

This result indicates that for large enough Reynolds numbers, the dissipation factor for 
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Figure 1: Solution of equation 3.3 



as a function of the amplitude o: 
of the classic Kolmogorov forcing). Dotted line: fx = 
fx = sin(27r?7) + A2 sin(2 x 2'kvi); Dash-dot-dot line: /; 
solid line: fx — sin(27r?7) + ^4 sin(4 x 27777). 



normalized by (3i = /36(oo)|/(x)=asin(27r,,),ae5R plotted 
the secondary mode forced, Ay^ (or a in the case 
a sin(27rr7), a G 5R; Dashed line: 
; = sin(27r?7) -|- A3 sin(3 x 27rr;); 



the Kolmogorov forcing is bound from above by a constant value equal to tt^/a/M. This 
value holds whatever is the amplitude of the Kolmogorov forcing. 

Another trivial case appears when the secondary term becomes dominant (i.e. A^'^ 1). 
The contribution of the primary term to the forcing function can therefore be neglected 
and the forcing becomes equivalent to a Kolmogorov forcing at a wave number greater 
than one: 

/3,(oo) = 4^ ^ ... : : — = (3.5) 



The upper bound on /? in the asymptotic case Re — > 00 is linearly proportional to the 
wave number of the largely dominant term in the forcing function. It is not surprising 
that when a single mode is largely dominant in the force shape functional, the bound on 
the dissipation factor /3 is related to this particular mode. However, the linear increase 
of /3b(cxD) as a function of the dominant wave number is not intuitive. This result implies 
that for a trivial force-shape, in other words a force-shape dictated by a single wave 
number, we could predict precisely the dissipation factor in the ideal case of an infinite 
Reynolds number. Each wave number forced alone is bound by a characteristic (3(, when 
Re — !■ 00. 

To further understand the force-shape dependence of the bound on the dissipation 
factor behavior, a complete resolution of equation |3.3| is required. Figure [l] shows the 
evolution of the bound on the dissipation factor in the asymptotic case, Ph{oo), for a non- 
trivial force-shape normalized by the value of /3b(oo) obtained for a classic Kolmogorov 
forcing {i.e., 7r'^/\/54), as a function of the amplitude of the secondary term. In addition to 



the properties exposed by equations 3.4 and|3.5| this figure also shows that /3b(oo) is very 



sensitive to the change in shape of the forcing function. Indeed, as soon as a secondary 
term is added to the original Kolmogorov forcing, /3b(oo) increases in a quasi-linear fashion 
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Figure 2: Dissipation factor, (3, obtained by DNS and normalized by /3i = 
/3|/(x)=Qsin(27r?)).QeSR plotted as a function of the amplitude of the secondary mode 
forced, Ak (or a in the case of the classic Kolmogorov forcing). Squares, dotted line: 
f = a sin(27rr/), a = 1,3,7; diamonds, dashed line: / = sin(27rry) + A2 sin(2 x 27r?7); 
down triangles, dash-dot-dot line: / = sin(27r?7) -I- A3 sin(3 x tttj); up triangle, solid line: 
/ = sin(27r7]) -|- A4 sin(4 x 2nri). 



as a function of the amplitude of the secondary term. This quasi-linear increase lasts until 
the rate of variation of /3f,(oo) is greater than about 80%. The increase of Pb{oo) then 
drastically slows down and plateaus. One can consider that the secondary term is largely 
dominant beyond this point. We can further anticipate that the characteristic features 
of the secondary term are the most obvious in the force profile, the contribution of the 
primary term being barely visible. We observe that (3b{oo) increases the fastest when 
A3sin(3 X 27r?7) is added to the classic Kolmogorov forcing, for small values of A3. We 
verify however that this feature is isolated and not linked to the odd wave numbers. 



4. Comparison to Direct Numerical Simulations 

Our Direct Numerical Simulations are performed using a fully de-aliased spectral code. 
The time stepping is based on a third-order Runge-Kutta algorithm for the non-linear 
and forcing terms. The viscous term is integrated through an analytic factor. The stability 
is ensured by a CFL condition. The flow is solved in a 27r in length periodic cubic box 
with 128 X 128 x 128 grid points. The viscosity is set to ly = 0.015625, large enough to 
ensure full spectrum resolution of the turbulent flow in our domain. 

Figure [2] shows the dissipation factor /? normalized by the average value of /? in the case 
of the classic Kolmogorov forcing at different amplitudes as a function of the amplitude 
of the secondary term. Each symbol represents a simulation with a different force shape 
following the definition of equation |3.1[ The symbols designing a particular secondary 
wave numder k are linked by straight lines using the same nomenclature as in figure 
[l]in order to ease comparisons between the behaviors of /3f,(oo) and /3. We verify that 
all the simulations satisfy the convergence criterion for the Kolmogorov flow defined 
by [Sarris et al. (20071. The simulations with secondary term amplitudes of Ak = 18 
verify kmaxV > where k^ax is the largest wave number of the simulation and r] 
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the Kolmogorov length scale, slightly under the commonly used flow resolution criterion 
kmaxV > 1-5- 

On figure |2] we recognize the same patterns displayed on figure [T] First, (H has a nearly 
constant value when only the wave number fc = 1 is forced. In this particular case, a 
limited number of simulations was possible as for an amplitude of 7 for the forcing, the 
limit of resolution of our fiow domain was already achieved. As regard to the force shape, 
we can anticipate that the value of P will remain about constant for any amplitude. 

Second, in the case of a non-trivial force-shape, (3 increases in a nearly linear fashion as 
a function of Ak until the secondary term becomes dominant (Ak sufliciently large). Then, 
the increase in f3 slows down and plateaus. The analytic bound on the dissipation factor in 
the asymptotic case of an infinite Reynolds number predicts qualitatively the behavior of 
the dissipation factor at our low Reynolds number Direct Numerical Simulations. Indeed, 
our simulations have Taylor-scale Reynolds number, Rx, varying between about 50 and 
100. The observations made from the predicted bound on the dissipation factor when 
Re — > oo apply at low to moderate Re: i) (3 is very sensitive to changes in the shape of 
the forcing function, ii) (3 saturates at a characteristic level for a given wave number and 
Hi) it seems that the level of saturation of (3 is linearly proportional to the wave number. 

The symmetry in the qualitative results obtain from figure [T] and [2] proves that (3 is not 
only sensitive to the force-shape as observed above, but /? strongly depends on it. First, 
when a single wave number is forced, a variation on the amplitude of the forcing function 
alone implies no variation on the shape {i.e., the general appearance of the profile) but 
a linear variation of the slopes of the force-shape as a function of the amplitude. The 
variations of amplitude affect the total kinetic energy of the flow, so it affects U (as 
defined in this paper, U cx E^^"^, where E is the total kinetic energy). The variations of 
the slopes affect the kinetic energy dissipation e. In addition, according to Kolmogorov's 
third law, e cx E^/"^ = e oc , for sufficiently high Reynolds numbers. Therefore, having 
/? (about) constant is directly related to the steadiness of the force-shape in this particular 
case. The same reasoning applies when a given wave number is largely dominant in the 
force-shape. Second, when a secondary term with a small enough amplitude is added to 
the forcing on a single wave number, the amplitude of the new forcing function varies 
barely whereas the force-shape becomes signiflcantly different (see case = 1 on figure 
[4c| for example) . The total kinetic energy of the fiow changes barely and the kinetic 
energy dissipation changes significantly. As a result, we observe a significant variation in 
P {(3(x e/£;3/2). We can anticipate the same type of infiuence on the dissipation factor 
for secondary terms with large amplitude, but not large enough for the secondary term 
to be largely dominant. 

As already observed, despite qualitative agreement, the numerical values for /? are 
about a factor 5 below the predicted upper bound (not visible here because of the nor- 
malization of /3 in the figures) (Doering et al. 2003 ). Also figures [l] and [2] clearly show that 
the theory does not predict correctly the magnitude of the increase in /? when forcing 
an additional wave number. The slopes in the quasi- linear part of the evolution of (3 axe 
greater in figure [l] than in figure |2] 



5. Mean velocity profile dependence on the force profile 

Figure [3] shows the mean velocity profile U{y), the Stokes fiow profile Ustokes{y), the 
optimal multiplier profile tpmiv) and the corresponding forcing function /(y) = sin{y) + 
sin{ky), where k = 2,3,4 respectively in a), b) and c). The Stokes fiow is the fiow 
associated with a lower bound on the dissipation factor (3. The Stokes fiow profile is 
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(c) 

Figure 3: Mean velocity profile U{y) (solid line and triangles), Stokes flow proflle 
Ustokesiy) (dashed line), optimal multiplier profile 5'„i(?y) (solid line) and corresponding 
forcing profile f{y) (dotted line), a) f{y) = sin{y) + sin{2y) , b) f{y) = sin{y) + sin{3y) , 
c) f{y) = siniy) + sin{iy). 



simply obtained by integrating twice the forcing function: 

Ustokesiy) = f([ F<i>{y")dy"\ dy' + C, (5.1) 
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(c) 

Figure 4: a) Optimal multiplier profile "tpmiv) for the force profile /(y) = sin{y) + 
Aksin{3y), Ak = 1,2,4,8. b) Mean velocity profile U{y) for the force profile f{y) = 
sin{y) + Aksin{3y), Ak = 1,2,4,8. c) Force profile f{y) = sin{y) + Aksin{3y), Ak = 
1,2,4,8, normalized by max(/), its maximum value. Triangles on dotted line: Ak = 1; 
squares on solid line: Ak = 2; diamonds on dotted line: Ak = 4; circles on solid line: 
Ak = S. 



where C is a constant always equal to in our configuration. The optimal multiplier 
in the case of infinite Reynolds numbers is given by: 

i^miy) = r sgn ( r <t>{y")dy" ) dy' + C, (5.2) 
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where C is a constant determined by using the zero mean condition on (j). We made here 
the hypothesis that the optimal multipher at finite Reynolds numbers is identical to the 
infinite Reynolds numbers optimal multiplier. 

In the cases presented on figure [3] the low wave number term in the forcing function 
is the dominant term. Contrary to the DNS presented by [Doering et al. ( [2003 j , our 
domain is periodic in all three directions. There are no free-slip boundary conditions 
to constrain the velocity profiles. Thus, the key feature of these velocity profiles is that 
they display all the characteristic elements of the force profile. These elements appear 
to various extent whether we have a Stokes flow or a fully developed turbulent flow. For 
the velocity profiles shown in these figures, the dominant shape is clearly a sine function 
because of the small contribution of the secondary term. But, we see by looking at the 
force profiles that a well defined change in slope induces a change in slope for the velocity 
profiles. These changes are more obvious in Stokes profiles (figure [s]) than in the turbulent 
flow mean velocity profiles. Indeed, in the turbulent case, a deviation in the force profile 
as seen on figure 3a around j/ = tt creates too little of a shear on a too little sub- 
domain compared to the combined effect of the shear on both sides of this sub-domain. 
Functions of high wave number provide less energetic contribution to the turbulent flow 
than functions of low wave number. Therefore, the weight on the secondary term of the 
forcing function must be greater than the weight on the primary term in order to have 
at least an equivalent contribution to the turbulent flow. Hence the minute contribution 
to the Stokes and turbulent mean velocity profiles of the secondary term in figure [3j The 
optimal multiplier profiles relate to the dominant term of the forcing functions. They 
appear to indicate only significant contributions to the force profile, here the primary 
term. 

Figure|4]allows a broader perspective on behavior of the optimal multiplier ipm - It shows 
the optimal multiplier profile, the mean velocity profile and the normalized force profile 
as the amplitude on the secondary mode is increased. As mentioned above, ipm appears 
to be less sensitive than the mean velocity profile as regard to the forcing profile. Indeed, 
we verify that for the cases f(y) — sin{y) and /(y) — sin(y) + sm(3y), ipm displays the 
same profile. For larger amplitudes {Ak ^ 2) on the secondary term, ipm displays a shape 
characteristic to the wave numbers composing the force-shape. As the amplitude on the 
secondary term increases, the amplitude oiipm decreases (figure 4a). The decrease in the 
amplitude of ipm is slower as the amplitude on the secondary term gets larger. We can 
relate the evolution of 4'm to the evolution of the slopes of the features characteristic to 
the secondary term in the force-shape profile. In figure 4c the local minimum at 7r/4 and 
local maximum at 37r/4 are caused by the addition of the secondary term. We observe 
that the variation of the slopes around these points is smaller as the amplitude on the 
secondary term increases. We verify that this variation of slope becomes constant as the 
secondary term becomes largely dominant (contribution of the primary term negligible) . 
Therefore, amplitude of the optimal multiplier remains constant when the secondary term 
is largely dominant. Finally, figure [4b| shows that the relation between ipm and the mean 
velocity profile is unclear. For f{y) = sin{y) + sin{Sy), the optimal multiplier shows no 
indication of the infiuence of the secondary term whereas the mean velocity profile does. 



6. Conclusions 

To summarize, low Reynolds number Direct Numerical Simulations confirmed quali- 
tative results predicted by the mathematical analysis at infinite Reynolds numbers. This 
confirms the dependence of the dissipation factor on the force-shape as predicted by the 
solution of the optimization problem. We can discern two major tendencies: i) when 
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the forcing function is mainly shaped by a single wave number, the dissipation factor is 
unique, its value being seemingly proportional to the dominant wave number, ii) (3 in- 
creases in a quasi-linear fashion when an additional term is added to the forcing function 
until this additional term becomes dominant. 

Also, we observed than the mean velocity profiles display all the characteristic features 
of the forcing function in a 3D-periodic domain. The extent of these features are modu- 
lated as a function of the contribution of each component on the forcing function. As a 
consequence, in SD-periodic domains, it should be possible to manipulate the force-shape 
in order to obtain a given mean velocity profile. The optimal multiplier evolution is tied 
to the force-shape, but its relation to the mean velocity profile remains unclear. 

In this short paper, we showed that in unbound turbulence, the force-shape plays a 
major role in the energy dissipation rate. The optimization problem solved by |Doering] 



et al. (20031 captures this dependence qualitatively for any forcing function. We demon- 
strated that their high Reynolds numbers prediction remains valid for moderate Reynolds 
numbers. Improvements in the quantitative prediction of the dissipation factor remain 
necessary in order to be able to use the dissipation in a turbulence model. 

The computational resources provided by the Vermont Advanced Computing Center, 
which is supported by NASA (Grant No. NNX 06AC88G), are gratefully acknowledged. 
We are also grateful to Professors Carati and Knaepen and all the Turbo team for pro- 
viding the code. 
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